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dynamic structure factors have been computed. A comparison with theoretical predictions was per- 
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I. INTRODUCTIONS 

A gas of interacting electrons imbedded in a uniform positive background is a basic theoretical model. Degeneracy 
(i.e. the influence of the Fermi statistics) and coupling (the role of interactions due to Coulomb forces) are important 
features of the electron plasma. The electron gas is characterized by 2 parameters - the density n and the temperature 
T (in energy units). The dimensionless parameter 6 = , Ep = ^(37r 2 n) 2 ^ 3 being the Fermi energy, describes the 
degeneracy of the electron system. The coupling constant of the electron plasma can be defined as the ratio of the 
average Coulomb energy to the average kinetic energy, 

7=^ = ^-, (1) 

2 Skin 0,T e ff 

where 

T eff = ^/ 2 TF 3/2 (/i/T) (2) 

is the effective temperature corresponding to the kinetic energy of the Fermi gas, [i is the chemical potential of the 
ideal electron gas, a = (3/47m) 1//3 and F 3 / 2 is the Fermi integral. 

For a strongly degenerate electron gas (# ^ 1) the coupling constant depends only on the density, 

Id = = 1.3575r s , (3) 

where rs = a/ clb is the Wigner-Seitz radius of the electron system in units of the Bohr radius. 
The coupling constant for an electron gas obeying classical statistics (8 ^> 1) is 

7c-r = ^ . (4) 

Important informations on the properties of the electron gas can be obtained from the knowledge of the dielec- 
tric function and the dynamic structure factor. The plasma dispersion relation, static correlation functions, and 
thermodynamic potentials can all be obtained from these functions. 

For a weakly coupled plasma (7 « 1, i.e. for very high electron densities or vice versa at very small densities 
and sufficient high temperature) the random phase approximation (RPA) jl]]^] gives good results in calculating the 
dielectric function because the kinetic energy is dominant. 

On the contrary, for electron densities corresponding to the conduction band density in real metals the average 
kinetic and potential energies of the electrons are of the same order of magnitude, and the plasma is so to say moderate 
coupled (7 ~ 1). The microscopic dynamics of electrons in such plasmas arc dominated by collisions and the non- 
collisional RPA becomes inapplicable. Exchange and correlation corrections to the RPA cannot be calculated exactly, 
but a number of approximate expressions of the dielectric function have been suggested H . 
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In order to check the validity of the different approaches microscopic simulations of the electron gas are of great 
interest. Classical simulations of the one component plasma were performed by Hansen et al. [Q. 

The aim of this paper is to report the results of a series of quasi-classical molecular dynamics simulations on 
the dynamic properties of the electron gas. The thermodynamic properties (mean energy) of the electron gas were 
investigated in a previous paper 0] . 



II. THE MODEL 



In order to treat the quantum electron gas by quasi-classical simulations we make use of effective pair potentials. 
Purely space-dependent effective potentials can be derived from the quantum-mechanical Slater sum |^,[7] . At short 
distances these potentials differ from the bare Coulomb potential and remain finite. On the basis of such potentials 
Norman and Valuev |8| and Hansen and Mc Donalds || performed molecular dynamics simulations of an electron- 
proton plasma. 

Another way to include quantum diffraction effects (i.e. the Heisenberg principle) into the dynamics is to blow up 
the phase space by introducing additional "quantum" degrees of freedom. This is done in the wave-packet dynamics 
developed to describe nuclear collisions |10|-|13|. However, this method leads to difficulties in describing thermal 
equilibrium properties of many particle systems . 

We use therefore in our calculations the ordinary 6N dimensional phase space, the particles interacting by effective 
pair potentials. However, a pseudopotential depending only on the space coordinates leads necessarily to the Maxwell 
momentum distribution. To model the momentum distribution of an electron gas governed by Fermi statistics we 
include in our simulations momentum-dependent interaction terms. We thus follows a line developed by a series of 
authors as e.g. Wilets and Kirschbaum, Dorso et. al. |lE|-|ls|]. 

In our simulations we substitute the quantum dynamics of the electron system by a phase space dynamics of 
Hamilton type with effective quasi classical Hamiltonian Q 

N 2 N N 

f-f 2m ^ r p ^ r p 

Here the first term is the ordinary (classical) kinetic energy of the electrons. The second contribution, the Pauli 
potential was chosen in a form suggested by Dorso et al. Jl5[ ], 

V p (p,r) = V e- A2 / 2 , (6) 

where A 2 = ^ + ^ is the effective phase space distance of two particles with relative momentum p and distance r. 

Po r o 

The last term in the effective Hamiltonian is the Coulomb interaction averaged with respect to the two particle 
Gaussian wave packets and is expressed by 

= erf(r/V2r Q ) 
r 

The Gaussian wave packet transforms the Heisenberg uncertainty condition into an identity, (Sp)(Sq) = h/2 which 
leads to 

r p Q = h (8) 

The other two parameters in the Hamiltonian Eq. (5) are chosen to describe the correlation function and the 
momentum distribution of a free Fermi gas. An appropriate choice is || 

h 2 

V Q = T eff , p 2 Q = mT eff , r 2 = — — , (9) 

ml eff 

with T e ff from Eq. (2), 

The simulations based on the Hamiltonian Eq. (5) with the parameters defined by Eq. @ result in a mean energy 
proved to be in good agreement with Quantum Monte Carlo simulations and with Pade approximations ||. 

That is the reason why we expand our considerations to the investigation of the dynamic properties of the electron 
gas within the developed approach. 

In restricting our calculations to a simple Hamiltonian given by Eq. (5) we make use from the fact that the collective 
dynamics of the electron system are dictated primiraly by the long range character of the Coulomb potential and are 
widely unaffected by the simplifications in the short range part of the effective potentials made in Eqs. (6) and (7). 
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III. DYNAMIC PROPERTIES 



Let 

N 

P (k,t) = exp(ik ■ fl(t)) (10) 

i=l 

be the Fourier component of the time-dependent microscopic electron density. The density-density dynamic structure 
factor is defined as the Fourier transform of the correlation function, 

i r°° 

S{k ' U})= 2^Nj 6 <P(M)P("M) > dt . (11) 

A closely related quality is the dielectric function e(k, lo) of the electron system. It describes the linear response 
of the plasma to an external electric field and is connected to the dynamic structure factor via the fluctuation - 
dissipation theorem (FDT), 

Qit. \ films~ l (k,uj) 

b{ ' W) = nirm [1 - exp(-/3M] ' 

where <f>(k) = 4ne 2 /k 2 , = 1/T . 

The imaginary part of the dielectric function is antysymmetric with respect to the frequency. Note that from Eq. 
(12) it follows that the dynamic structure factor does not possess symmetry, but satisfies the relation 

S(k, -lu) = e~ fihu S(k, u) . (13) 

The dynamic structure factor defined by Eqs. (11) and (10) is directly measurable in the MD simulations if one 
identifies the Heisenberg operator r|(i) with the position of the i-th particle in our simulations. 

However, the thus obtained quantity (we denote it by R(k,ui)) possesses symmetry. It corresponds therefore to a 
classical FDT, 

R(k,uj) = (mr^Piu)- 1 Ime~ 1 (k,uj) . (14) 

From Eq. (Q) one concludes that R(k, u>) can be regarded as a normalized loss function. 

The quantum-statistical dynamic structure factor obeying the relation Eq.(13) has to be calculated as 

S &") = - h f R(k,w) . (15) 

1 — exp(—pnuj) 



In what follows we will regard the normalized loss function. Note, that in the classical case the loss function and 
the dynamic structure factor coincide. 

In our molecular dynamic simulations we used the algorithm of Verlet |19| to integrate numerically the equations 
of motions obtained from the effective Hamiltonian of a system of 256 electrons. The typical length of the MD runs 
were about lO 3 ^" 1 (to p being the plasma frequency). 

The equilibrization phase was replaced by a Monte Carlo Simulation using the algorithm of Metropolis J2(| . 

The forces were calculated by an Ewald method in order to account for the long range of the Coulomb interaction 

13- 

The motion of the electrons can be studied by calculating the velocity autocorrelation function < v(t + r)v(t) > t . 

We see (Fig.l) that for T = 1 the velocity autocorrelation falls monotonically to zero , whereas for T = 100 the decay 
of the velocity ACF is characterized by appearance of oscillations with a frequency close to the plasma frequency. 
That means that the motion of a single electron is coupled to the collective density fluctuations. 
The collective motion is described by the dynamic structure factor (or the loss function). R(q,uj) (q = ka) is plotted 
for two q values, at T = 1 for different 9 = 1 (moderate degenerate) and 9 — 50 (classical) and for V = 100 at 9 = 50 
(strongly coupled, classical electron gas)(Figs.2-4). The results of the simulations are compared with theoretical 
predictions from RPA. At moderate coupling constants (r = 1) the plasmon peak of the loss function is less than that 
predicted by RPA and slightly shifted to the left (Figs. 2,3). In both cases the plasmon peak can be observed only for 
the smallest q value. The change of 9 in the range from 50 to 1 has only small influence on the results (Fig. 5). 
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In the strong coupling regime r = 100 the observed plasmon peak at the smallest q value is extremly sharp and 
centered close to up (Fig. 4). A well defined collective plasmon mode has been developed. At a q value 3 times larger 
the plasmon peak widens, but is still present and shifted to the left by about 10 per cent. At still larger q values the 
plasmon peak vanishes. 

This behavior is in striking contradiction to the RPA predictions where no plasmon peak can be observed due to 
the strong Landau damping. However, the RPA is inapplicable in the strong coupling regime, where the potential 
energy is dominant. On the contrary, our MD simulations for the case of weak degeneracy are in a good agreement 
with corresponding MD simulations of Hansen et al. for the classical one component plasma [Q. 

We conclude therefore that our model yields reasonable results in describing the dynamic properties of the electron 
gas at least at weak and moderate degeneracy. 

IV. ANALYSIS OF THE RESULTS 

Important characteristics are the frequency moments of the imaginary part of the inverse dielectric function(DF). 
They are defined by 

C v {k) = / uj"- 1 Ime- 1 (k,uj)du , u = 0,l,... . (16) 

J-oo 

Due to the antisymmetry of the imaginary part of the inverse DF all even frequency moments vanish, whereas 
the odd frequency moments are purely expressable in terms of the static properties of the electron gas. After a 
straightforward calculation one obtains, 

C (fc) = (l-e- 1 (fc,0)) , (17) 
C 2 (k)=u 2 p , (18) 
C 4 (k) = u*(l + K(k) + L(k)) , (19) 

where 



K(k) = 3{k/k D ) 2 + ^7r/18{\^k 2 /\ L ) + \ 2 T k A k D 2 (20) 

is the kinetic contribution involving quantum corrections, k 2 ^ = 47me 2 /3 , At = (h{3/2m)~ 1 / 2 , Xl = 3/2e 2 /3 and 
u 2 = Anne 2 /m . The contribution 

1 f°° 

m= 3^n~ I P 2 i S (p)~ 1 U(P,k)dp (21) 



o 



takes into account the electronic correlations, 



A quantitative analysis of the results of the MD simulations should be based on the frequency moments defined in 

The Nevanlinna formula of the classical theory of moments constructs a dielectric function which satisfies the known 
sum rules Co to C4 |f22f : 

e-\k,z) = 1+ . _ "f* + ® 2 2T ■ (23) 

z(z 2 -uj 2 ) + Q(z 2 v ; 

Here Q = Q(k,z) is an arbitrary function being analytic in the upper complex half-plane Imz > {lo = Rez) and 
possessing there a positive imaginary part, it also should satisfy the limiting condition: (q(k, z)j z) — > , as z — > 00 
within the sector § < arg(z) < ir — i3 (0 < $ < ir) . 

The frequencies uJi(k) and u)2{k) are defined as respective ratios of the moments C n {k): 
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^ = C 2 /C = ^(l- e -i(fc,0))- 1 , (24) 



u% = C A /C 2 = + K{k) + L{k)) . (25) 

There is no phenomenological basis for the choice of an unique Q(k, z) , which would provide an exact expression 
for £ _1 (fc, u) . 

If one is interested in the investigation of the dispersion relation only it sufficies to neglect Q(k, u)), since the 
damping is small in strongly coupled plasmas. If one puts Q(k,uj) = one obtains the expression of the inverse 
dielectric function obtained within the quasilocalized charge approach of Kalman j23|. The disadvantage of this 
choice of Q is that damping is not taken into account. As a result the shape of the dynamic structure factor within 
this approach is reduced to a simple 5 function peak at the frequency u>2{k) . Thus only the peak position but not 
the shape of the dynamic structure factor can be described appropriately. 

The easiest way to go beyond the simple approximation Q = is to put the function Q(k,oS) equal to its static 
value: 

Q(k,z) =Q(k,0) =ih(k) , (26) 
where h(k) is connected to the static value of the dynamic structure factor S(k, 0): 

h{k) = w [ (W2/Wl)2 - ^ ■ (27) 

k d S{k,0) 



From the Nevanlinna formula the loss function reads 

R{k,u) = S(k,0)- r ^iM , 2Nj 

Formula Eq. ( p8| ) interpolates between the low frequency behavior (described by Co(fc) and S(k,0)) and the high 
frequency behavior (given by the moments C^fc) and C^(k)) of the loss function. Though there is no direct justification 
that the above interpolation formula describes the loss function in an appropriate manner also for intermediate 
frequencies, we expect that Eq. (^8|) reproduces the whole shape of the loss function at least qualitatively. 

The analysis of the MD calculations is based on the formula Eq. (p8|). To calculate the frequencies h(q) , u>i(q) and 
LU2(q) (q = ka) the static structure factors from HNC equations were used. Since the latter are classical equations we 
have restricted our comparison of the sum rule approach to the MD calculations to the case of weak degeneracy. (In 
this case the dynamic structure factor S(q,u) coincides with the loss function R(q,iu)). 

The results of the comparison of the loss function calculated from the MD simulations and that obtained from the 
Nevanlinna formula are shown for different q vectors at T — 100 and r = 1, respectively. The agreement between 
theory and simulations is quite good. The theoretical curves reproduce rather well the variation of the shape of the 
dynamic structure factor and describe the plasmon peak position in a good manner. However, the agreement of the 
height of the peaks is less satisfactory (Figs. 6-8). We believe that one of the reasons for this disagreement between the 
results of simulations with theoretical predictions might be the normalization to S(q, 0) which is a value rather bad 
measured in the simulations due to the poor statistics at long times. 



V. CONCLUSIONS 



The results of quasiclassical molecular dynamics simulations of the electron gas using momentum dependent effective 
potentials have been reported. The quasiclassical MD computations were performed for different coupling constants 
(F = 1 and r = 100) at various degeneracy (9 = 1 and 9 = 50). The effective potential was chosen to describe both 
the Pauli principle and the Heisenberg uncertainty. The changing of the coupling constant T leads to qualitative 
changes in the dispersion curve of the density fluctuations. 

At moderate coupling (r = 1) the dispersion is positive, the plasmon peak is observed only for the smallest 
wavenumber q = 0.618 . As can be seen from Figs. 2 and 3 this behavior is qualitatively confirmed by the RPA 
calculations, however, quantitative deviations from the RPA predictions as a slight shift of the plasmon peak to the 
left were observed. 
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In the strong coupling regime the shape of the loss function changes qualitatively (Fig. 4). At T — 100 the dispersion 
is negative, a very sharp plasmon peak is observed for the smallest wavcnumber, with increasing q the peak widens 
but is present up to q = 3.1. At still greater q the plasmon peak vanishes. 

Thus, the collective behaviour of the electron gas at weak and moderate coupling can be understood as oscillations 
of the total charge in the Dcbye sphere, the individual particles in the sphere moving almost independently. This is 
the regime of the collisionless plasma described by the RPA. On the other hand one can interpretate the collective 
motion in the strong coupling regime as solid-like collective oscillations when the motion of each of the particles is 
coupled to the collective oscillations. This is also confirmed by the shape of the velocity autocorrelation function 
(Fig.l). 

In contrast to the qualitative change of the shape of the loss function by varying the coupling constant the variation 
of 9 in the range from 9 = 50 to 9 = 1 has only small influence on the behavior of the loss function (Fig. 5). However, 
we expect a greater influence at higher degrees of degeneracy. 

As demonstrated in Figs. 6-8 the above features are rather well produced by a simple sum rule analysis. 

Finally we note that our quantum molecular dynamic simulations describe the dynamic properties of the electron 
gas only approximately. 
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FIG. 1. velocity autocorrelation function for 9 = 50 at different!" (1-F = 100, 2-T = 1) 
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FIG. 2. comparison of the MD loss function R(q,co) versus frequency uj/lo p with the corresponding loss function from the 
RPA for different wavevectors q at Y = 1 and 8 = 1 
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FIG. 3. same as Fig.2; for T = 1 and 9 = 50 
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FIG. 4. The MD loss function R(q,ui) versus frequency lu/lo p for different wavevectors q at T = 100 and 9 — 50 



8.0 



6.0 



o o thcta-1 

theta=50 ? 



cd 
DC 
i 
DC 



4.0 



2.0 



P : ■ 



o--o" 



0.0 



0.0 



0.5 1.0 1.5 

w ( in units of the plasma frequency ) 



2.0 



FIG. 5. The MD loss function R(q,uj) versus frequency lu/u p for wavevector q = 0.619 at fixed T = 1 and different 
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FIG. 6. comparison of the MD loss function R{q,u>) versus frequency ui/lo p with the corresponding loss function from the 
sum rules approach (Eq.29) at T = 100 and 6 = 50 for wavevector q = 1.856 
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FIG. 7. same as Fig.6; at V = 100 and 9 = 50 for wavevector g = 3.094 
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FIG. 8. same as Fig.6; at T = 1 and 6 = 50 for wavevector q = 0.619 
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List of figure captions 

Fig.l. Velocity autocorrelation function for 9 = 50 at different T (l-T = 100, 2-r = 1). 

Fig. 2. Comparison of the MD loss function R(q,uj) versus frequency lo/uj p with the corresponding loss function 
from the RPA for different wavevectors q at Y = 1 and = 1. 

Fig.3. Same as Fig.2; for T = 1 and 9 = 50. 

Fig. 4. The MD loss function R(q,Ld) versus frequency uj/lo p for different wavevectors q at T = 100 and 9 = 50. 

Fig. 5. The MD loss function R{q 1 u) versus frequency uj/oj p for wavevector q = 0.619 at fixed T = 1 and differ- 
ent 9. 

Fig. 6. Comparison of the MD loss function R(q 7 uj) versus frequency lj/uj p with the corresponding loss function 
from the sum rules approach (Eq.29) at T = 100 and 9 = 50 for wavevector q — 1.856. 

Fig.7. Same as Fig.6; at T — 100 and 8 = 50 for wavevector q — 3.094. 

Fig. 8. Same as Fig.6; at T = 1 and 9 = 50 for wavevector q = 0.619. 
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